{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Beam Expander Design"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "from optiland import optic, optimization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define a starting lens:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0QAAAChCAYAAADneH+CAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/TGe4hAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAw50lEQVR4nO3dfZQU5Z0v8G+9dFX3vCJvM46AYmRXUaNEIkv03JgNV9Z1YxL3mMRDyET3xCOLJxq8Rt2Nupo1qLkxRtejiWeN2ZtEE89qVt1jDEGicUVAkBhF0awsqDAg4DAzPdPd9fLcP7qquqq6uqdnpnuqe/r7OYdD11vPw5Qw/fX3e56ShBACRERERERETUiOewBERERERERxYSAiIiIiIqKmxUBERERERERNi4GIiIiIiIiaFgMRERERERE1LQYiIiIiIiJqWgxERERERETUtNS4BzBRtm1j7969aG9vhyRJcQ+HiIiIiIhiIoTA4OAgenp6IMuV1X4aPhDt3bsXc+fOjXsYRERERERUJ959913MmTOnonMbPhC1t7cDAP77nV3eayIiIiIiaj6Dg4P4yPHzx5QLGj4QuW1y7e3t6OjoiHk0REREREQUt7FMpeGiCkRERERE1LQavkJERERERETxEELAMAwAQCKRaMhFzhiIiIiIiIimKNsWyJo2MqaFjGEjY1j5bSO/nTUtjHivneOGnd9n5l9nTBtZw0LGtJHJWd57ZU0LuZyBv8y+CABYc/X/gaZpMf+Jx46BiIiIiIhokgjhBBQnnHjhwg0cTjgphJHC8REnhPiPjxiF8/znZ5z3y5l2xWOTJCCpytATCpKqjGRCyf9SZegJGUlVQauuYEabhmRCga7K0GUB85UafsMmQU0D0fPPP4/vfve72Lp1K/bt24fHH38cn/vc57zjQgjcdNNNeOCBB9Df34+zzjoL9913HxYsWFDLYRERERERAch/Hs2ZdnQYCYWUbKjK4g8jXrUl56vAmLbveKEKMxZJJ4i4gSSZkL0wkkwoaNEUHNWiIZUIBhldlZHSCucFgo4WPC+ZUJBKKEgo0phb3nK5HO5kICotnU7jtNNOw6WXXooLL7yw6Pgdd9yBu+++Gz/5yU8wf/583HDDDVi+fDl27NiBZDJZy6ERERERUR0SQsCwRCiMWIW2rVBQ8bd/Bc7znT/iOy8bDj6mDSEqH18+QMjQfeEkHC46WxJIJRRfkMmHkVRCCVzrDyP+wKOrClIJGZoqN+ScnEZT00B03nnn4bzzzos8JoTAXXfdhW9961v47Gc/CwD4t3/7N3R1deFXv/oVvvSlL9VyaERERERUATegZM2IMBLRopX1hY2iOSum2/pVmJsSrLbkzxlLQNFUOdTm5YYUJ2QkFLQnVaf9yxdknOOlwoieUIqqM7oqQ5YZUKaa2OYQ7dq1C319fVi2bJm3r7OzE0uWLMHGjRtLBqJsNotsNuttDwwM1HysRERERPXEsAqtWEXhItSiVWj/Ck+gLz83xR9o7DEElIQiFbdohcLFzDa9aG6KF0a8Nq9gaPFXYPzHGVBoomILRH19fQCArq6uwP6uri7vWJS1a9fi5ptvrunYiIiIiMbCtEIrcRkR1ZHQBPlyE+jDc1bc393WL2sMCSWhSEUtWoE2r4SM6a1aoP2rKMgkiueiRM1ZSSYUKAwo1GAabpW566+/HmvWrPG2BwYGMHfu3BhHRERERPXGskUgjIyEKiKZUItW8QR6q7ASWGA1LyuizcuGOYaAosiSV90oatFywsa0lkTkBHr/hPjoCfRuJaZwXFXkGn6niRpfbIGou7sbALB//34cffTR3v79+/fj9NNPL3mdruvQdb3WwyMiImp6Q0NDGBkZqcp7mXZ+Ja+sZSNr5l9rqVaYQh61OhI1gb7wnJTgBHo3+BhW5QFFllDU0pVKBMNFe1LFrHbNm5sSmBjvtnT555p4bV7Fc1YSDChEdSW2QDR//nx0d3dj/fr1XgAaGBjApk2bsGrVqriGRURE1FQsW2A4ZyKdtZDOmRjKmEjnLHzQP4innn4GhpBgQYYlZFiQYUKG5ewzvf2S73XxtgkJApWHADeg6KFWrJQvXLQmVe9ZKKWWJQ6HkUJlJVidGc9Sw0Q0ddQ0EA0NDeFPf/qTt71r1y5s374d06dPx7x583DVVVfhn//5n7FgwQJv2e2enp7As4qIiIgoKGfaSDshZihrIp3Nh5h01nS2reBxN/CEzk3nLAznrDJfaT4kAJoCaIrk/UrI+e2k77V3zD1XlqAp+fkr4es0RcKut3dCFhYu/uJFRe1fDChENJlqGohefvllfOpTn/K23bk/vb29eOihh/DNb34T6XQal112Gfr7+3H22Wfj17/+NZ9BREREU4oQAhkjH2K8wOIEEi/QuBWaqEATCjyjtYOlEjJadRVtuopWTUGr8/vsDh2tWitadQWtmnNcLxzPb+f3pfsP498f/n/4m/PPw8yZM6v+PdmW7cN7772HBbPbqv7eRERjUdNAdM4550CUWUhekiTccsstuOWWW2o5DCIiojGzbYHhXHGFZchfYQlXYHwtZ/7AM5yzyq4KJklwAkohnLiBZnprCq2a6mwXH/cCj3O8RVOrsspXX3YAiiRqVqlRFAWWVa46RUQ0ORpulTkiIqJS/K1kY2opC4ScSlrJ8q1gUeGkI5XA0dNS+fCiFSowUeHFPZ5KKE3XIiYrCkzTjHsYREQMREREFB9/K1nJCksmVIEp1XKWs5Az7bJfz20lC1djZrfrXqtYVCuZF2h8xzWVK4VNhKIosO3y94uIaDIwEBER0Zi4rWT5cBKsqoTnw0RN6g+0nFXYStaqK0UVFn8rWaumoC3pny/jCzzO8RZN4fNY6ogiy6wQEVFdYCAiImoChmWXqLD4Q03UpH7fcSfEVNRKpvkrLPnX7UkV3Z1JX8UleLzNF17c483YStYsZKdCJETt5ikREVWCgYiIqA4JIZA17aJwUqi4FCozQ75qTHEFJn98tFayZEJ2KjD+tjAFs9p0HDcjum0sam4MW8moUrKc/+/EsiyoKj+OEFF8+C8QEVGVFLWS5UykM75WspyJdMYqMx8mWI2ptJUsXGGZ05KKmMSvBAOP7zhbySgOCgMREdUJ/gtERE3N30oWbA/zr0AWOl6iWjNaK5kqS6EKixNUkiq6OvSIlrHowNOq5VvJ5CosrUwUF1lRAIBLbxNR7BiIiKihuK1kpR9qGTXJP6oCk9+fraCVLLDqmBNOZoZayQqT+n0tZ8ngcU2VOVeCyOFvmSMiihMDERHVnG0LDBtWmQqLL7SUO+4EmtFayVo0pXg+TFLFMU4rWavmq7aUWGa5TVfZSkZUQ24g4tLbRBQ3BiIiimRadokKzNgfejliWBClM0xEK5kTTJIqZnfoY5oPw1YyosagsEJERHWCgYhoihBCIOeuSpazSlZgws+AKfXQy0pbycIVlpltOo6dHrWEcun5MGwlI2o+EgMREdUJBiKiGPlbycpVWNKBkFP6oZdmBa1kwYdc5l/3TEuVr8BEzIdhKxkRTYTCljkiqhMMRERj5LaSjVZhKT3pv/B6ODd6K1nRHBcnpLitZJVUYFp1FS1sJSOiOuJWiGxWiIgoZgxENOW5rWRlw0m4QlOq5SxnImOU/7+ZuipHVlhmtGmYp41tUj9byYhoqvJWmWOFiKghCSFg2gIZo/H/pwYDEdUlIfIPuBzrfJh0iWpNuVYywLcqWajCcvS0ZPkKTGhSf4umIMFWMiKiURVWmWv8D1NEUdzAYFoCli1g2DYsZ9v/2rTzvyzLhuE737Ttouuj3qvwHvl9hu96o8R7mc7XCo6h3Pjs4FidXwCgwsLKVMzf7AliICLPjh2v49DBQxWdawsBwwZyFpzf89vh18XHgLZp0yFUvfR8mQpayRRZQpt/RTJfhWV2hx5cqSyyAqPmr2crGRFRLNzqt22xQtSs3MBQ+oN34QO+PzAUzi8ODCXfyxKw3PPLXO8PDP4P/4YVDAzF4wuHEVH2ERETpcoSVEWCIktIyLL3WpUlqL7thCxBVWTfscJ2iyZDlfOvExHXq7LzHhHXu8dVWYJkW3jzN6/U7M86GRiIqmjdjv1Y98aBuIcxLoZhYOfOnbCEBKEkYEGGKSRYQoYJKfRaho2xBQgZNlRJQBEWdDWD7plHeYFmequGudPDk/rVYOAJLbuss5WMiKih8TlEoxOi8OG6oiqAFf2Bv2zgsIQvZASvH7WqYdmBwBAVOMJj9f9ZJjMwKM52VGBQ3A/5ocCQ0twgIEcGBnc7KjC4ISMfKoIBwg0ZRYFDcd47MNbC9e5YFVmqq89AuVwOb/4m7lFMDANRFQ1mTeztz8Q9jHGxbAtDQsOs6dPQ2dYCTZGgKRJ0VYauStAUGZoqQVfC2/lzNEWCpsqh7fz5mipBdv7ibnzxRaSHh9Hb+9cx/4mJiChO1ZhDFAgMoSpAqVYhI/Ah3i66vuhDfKAqEawalGot8qoGRe/vHCvRFlXU1lT3gUGJITAEKxilrq+nwED1j4Goii5cdAwuXHRM3MMYFyEE7rj9Nnz8lDPx538+N+7hEBFRDbgtSjkzX1XIWbbz2kbOcvaZ+X0HDg5ij9WJTe9n0NJ/BKaN/LwBr80JXjCwReF1/oM8YAl3u/Da9G3nDAvvZRdg27p90H/fXxQYKm2LqpVqBwa1xPWBNqSIwKBEBoBgyIgKDKMFDgYGogIGIgKQ7+VuaWlBZmQk7qEQETU8MxQw3MDhBZGofZazz3SuNSP2ha71B5iifb5r/e9fbn5msQVYv/EIgCOBvbIE74O14n8tw9nnfvBG8LV3jQRJAhTYaNMkTO9IjjInwhcSSgUGXyAYfd4DAwMRFTAQkae1tRWZTGO2/BFR83FbngphIBQWTDdwCF9YcMNF6Nqoff5rKwkwvq8/kcKFLAGaKiOhyNCUfCuSprqvZeeY27osQ1dltCfVwL7Ia9XK3s/d13/oIP790Uew/Nz/jVkzpnvhQpZQlbBgGAYeefh5XPCJU7Bw4ckTfj8iovFiICJPCwMREYXYTugYrb1q1OpEiYpFcYAJV0DyoSR4bf68iU7I1sYYEFp0pXzg8F9bwfuVCjBqnSzdr+U0pCQTbZqMFk2p+vt7q8zVsO2NiKgSDETkaUm14NChg3EPg6jpCCHKBw5ftWNCgcMNG/5zLRHcDl070TkaCUXyBYDCAiyFAJBfeMWtYLQnFO98reja/Lkl9/kCSLnAoan59iq2RMXLXVRBcJU5IooZAxF5Ui0pZPdm4x4GUU1ETSYvOTfDDFYngq1SJfb52q0qCjC+aw1rYqFDlaVRA4C3rcpo1RVMUxPlKxiKjIRaaQuWL9R4YSZ/LkMHjcYe26QmIqKqYyAiTyqZQjaXi3sY1OBMayyBYwITzsfyfuOaTB7kTsCO+tAfrFbk96U0BR1KIrCv1LUTacFKKDIfLEwNyQ3LrBARUdwYiMiTTCWRy7JC1AhsOxwQoiaIjxIQzOoHjpw5scnkkoSSFYeoditNldGmq9EtU+F2qzG1YBUHGIWhg6iqJClfQRSsEBFRzBiIyJPUk7AsC5ZlQVGqP4G20dh2ocWqVPgo2xplOcvmRu0LTxAfbYWr0Nef6GTyUlWNUu1RnS2JiU8i91qsSs8h4VK3RM1FkiTYghUiIooXAxF59KQOAMjlckilUpPyNf2TyUtOBq/gAYKVPH8j6v3KTS6f6LwObzK5LwgEt4PtVm1JpWwLVqkJ4yXbt0q0dHFeBxHVC1aIiKgeMBBV0eZdh7Fl94dxD2Pc+vsH8UfjaKS3HoSa0LwnjZu2gO17bfmfUC6cp457r0XoaebBYyOZDph2O35w63ovmEyEIkvlV60KtVu1agqmpRKRE8cDASRqnxtAVDlU+SgOMKoscV4HEdEoJEmC4LLbRBQzBqIqem3vAB7e/G7cwxg3y7YxYs7CnnfS0NSM74ni8B7Ip0Q8dVxTJCgJGWroCeVyxLn73nsXpmngLz5+0qiT0kuFD3/I4bwOIqLGJUkSBBiIiCheDERVdOlZx+HSs46LexjjdujQITzwox/i3OXL0dXVVZOvsdH4E9LDw+ht4O8TERFVD1vmiChu9fE4bKoLiUQCAGCaZswjISKiZsA5RERUDxiIyMNAREREk0mSJLBjjojixkBEHlXNd1BalhXzSIiIqBmwQkRE9YCBiDwMRERENNm4qAIRxY2BiDySJEGWZdgMRERENFmYh4goZgxEFCDLMmybTw0nIqLa47LbRFQP6iIQ3XvvvTjuuOOQTCaxZMkSbN68Oe4hNS1FUWCzn5uIiCYLf+QQUcxiD0S/+MUvsGbNGtx0003Ytm0bTjvtNCxfvhwHDhyIe2hNiRNciYhoMrFCRERxiz0Q3Xnnnfja176GSy65BAsXLsT999+PlpYWPPjgg3EPrSkxEBER0aSRpLhHQEQENc4vnsvlsHXrVlx//fXePlmWsWzZMmzcuDHymmw2i2w2620PDAzUfJyVGslZeOdgOu5hTMhBK4XUoID+wUhN3v/9EQWZbAKv762f+0ZERMUOHRrGQbsF73xo4JBdm58JH5hJ7B6wcRR/JhA1LNMwvNcZw4KmxTiYcYo1EB08eBCWZaGrqyuwv6urC2+++WbkNWvXrsXNN988GcMbs3cOpnHh/S/FPYwJOh74I4A/7qrR+7cDaMeDDf99IiJqBgvx5PrDAA7X6P2PB17KAi/xZwJRo1JhYWUq/3rXwWGc1pqKd0DjEGsgGo/rr78ea9as8bYHBgYwd+7cGEdUcPzMVjx2+V/EPYwJeeThn2PevHk4YcGCmrz/q6++ikwmg7/5m8/U5P2JiKg6Dh06iCeeeAJnn302Ojs7a/I11q9fj5NOOgmnn76oJu9PRLVnGgae/vkrAID5M1tiHs34xBqIZs6cCUVRsH///sD+/fv3o7u7O/IaXdeh6/pkDG/MUpqCk3s64h7GhMxURjC3XcKfzapNuj+UspAWRsN/n4iIpro+eRgz5WEcf1QCM2bU5mfCq2oGx3bI/JlA1MByuRyedl4nE0qsYxmvWBdV0DQNZ5xxBtavX+/ts20b69evx9KlS2McWfOybQGJk1yJiGgyCAEJ/JlDRPGKvWVuzZo16O3txeLFi3HmmWfirrvuQjqdxiWXXBL30JqSEDZkBiIiIiIiahKxB6IvfvGL+OCDD3DjjTeir68Pp59+On79618XLbRAk8OyLEhy7KuxExFRExBCgAUiIopb7IEIAK644gpcccUVcQ+DkG9ZVBiIiIhokrBljojixk++5LFtG0IIyEpjTogjIqIGxDxERDFjICKPaZoAAIWBiIiIJoEQXMiHiOLHQEQe08w/aZiBiIiIJgtb5ogobgxE5DGMfIVIVetiahkREU1xrBARUT3gJ98qOjCYxcGhbNzDGLf+D/txyE5hbxowPsxCkQFVlqBIEhRZgiIDiiTl9znbXKKbiIjGi6vMEVE9YCCqol++/B7u2fDfcQ9jgk7GE8/2A+iv6GxZghOY4IQkCYqEQmiSJKjOMVmWkEm3ASKF9T/eAk2RoakyEv7fFcn32j1WYp8SulbNXx98v/y5CUWGKkv8P5FERHWG/y4TUdwYiKroC4vn4FN/PivuYYzb3vffxzO/eQb/65PnQE+1wLIFLFvAtAFLFL/O/8pvm0WvnW1bOOcDpi2wd98RZE0Ls9p05Cwbhmlj0DCRM20Ylp3fZ9nImcL32j2Wf9/xkiQUhSR/qCoKU6p7XAoFLGef/9qofeFry7yfIvMDARE1H7bMEVE9YCCqotntOma363EPY9zUwb2YKQ/jlKPboCeTNfkaG1/chfTwMHov+ui4rrdsURSS/GHKe22K4n1W4drIfVboWtNGOmv6vp4IvZ+zz7lWjD+rQZYweugaQ8DyrlUn/n4ywxoR1QgDERHVAwYi8mQy+flPCU2LeSSl5dvyFCQT9bUSnnCqZuGQFAhTUfvCAcssH7rca7OmjcGMOfr7OddOhCpLlVXJIgKWV4VTI86N2ucPct6+4lZK9+vzgxRRY7NtG7LE9Z2IKF4MROTJZjNIJBKQZf5wGitJkqAqElSl/r53QggYligbsCoKbGbw+mAFrtDimM5ZyJlGYF+gIuf7OoY1gbIagITia30sMYfMDWmR+wLz0KSIeW2hfRW+H8MaUeX4d4WI4sZARJ7MyAi0Oq4O0fhIkuQtTNEa92BCbDs/5ywcxnJWqOXRH8YqqapFzUEzbaQNEx9WUIEzrPy4JqKi0DXWwFbJ+6llFidR8vPV+AGU6oEQIt8yx7ZcIooZAxF5RkZGoOuNOweKGo8sS9Dk/If3euPOVytXMRttXplhlg5d4WrbkWGjogpczrQxkawmSSi/WmM4THnz0CpdiKT0So+jtV1ycZHmIpyJlxJb5ogoZgxE5BlmICLy1Ot8NQAwreKQVBSmRmtxHCV0+a9NZ82KWiYnOl9NceerlZubFqqYVXUhkjLvx8VFqs8NRGzTJqK4MRCRZ3h4mIGIqAGoigxVAVKor7AmRKEFsvLAVmJVyNFWirRsjOQsHLGMiubATXS+mipLlYUu/3y2Ei2OFS0uUhTiQhW9KbC4iBeIGnT8RDR1MBCRZzidxqxZjfscJSKKlyRJ3kIX9cZdXCTc4ugPXm6LY8lAZ44S8rz3EEhnLeSGjaKAFnVtNearjba4SKVti/5ANpIewlvmDLTuHkHnh0cCD9oOvPY9lFuRJefB3M4DuiX3OLzX7nFWiIioXjAQkSedTmPuvHlxD4OIqOr8i4ugzgrhtl28CuSYWxwrqco5+4dCD8MuHQBtWGI+/mvLAICBqv+5JQASPoafP/o+tERfvmXSCVeqIkN1wpPqLAaSkPOreeZfy855hXNU59qEc77/+sIxN5TJRdf6r1e9cbjvI0NR3DHIZcbqvLcssc2SqIEwEFVR/v/0TayHPi6WZWNwJAdNq7NPCkREU5wsS9BlBXqdzVfr6+vDgw8+iOXn/TU6j5oOy84/b80SgOm8Np1t71jouGU7r0X+teW9FhjJGtiydStO/ejpmD5jFkw7H8zc9zVtG6Zv27Bs3zEB08ofzxo2DNvIH7NCx90xhrb9ryfyUO1yJAlFocwf8FQvoAUDlRv43NejBzbZF8ac8OYPfyUCm3t9VLgsBNBwuGTgo6mJgaiKfvT7Xbhnw3/HPYwJ+Bh++jsTynM7nJYLKb/ilALogW3nlwzvddExBYHthCJBVyTsPDCMZELB/oEMWnUVLQmF/6gSEdUpWcq35KUS1W9rGx4eRvoPH+DzJx+FE044rurvXyk7FMD8ocmyBYzQa8tyzvUFMKtUYPPOca73B75QYDOcbSt0bdYwYfjCoVEi4JU6VqvAJ0sorsgVhbdC4AtU53zhT1XkktW2yGpgqLpXKrBFXV9x9ZCBr+kwEFXRuQtnY/7MlriHMS65bA7PPPNM/h9TyDCFDMuSYFkyTMiwhAwLEoYhY0DIsCA750n5351ty7edb4gI6wEA/Oj/Pu/tadEUtOoq2nQFrZqKVt3dVtHqHGvVVLQlfcc157heON6qK3U5d4GIiIoJp6Mi7jlE3vL/mJo/P6yIwBYV3grhT8AKhcNwYAtcbxWHSTeURQU0/zHTEsgYZqAa6A+gURU9b3xOwKwVWUJRwCuutkVU5CICW6Hds3w1sHwraHGgLAp/JaqH7vX+sTbqYiy1wkBURX/W1Y4/62qPexjj9pcLvoj00FBV3stdbSrrTELO/y6QNW3obR3I2RLSWQtDWRPprIl0zkI6a+a3ndf7+jPOtol01kI6ZyJjlG9J1FW5ODBpaiFw+UJW8LjibBeO66rMfzCIiGrE9hZV4L+zteQuelGPz3urhqjAFQx/hQAVCFRW6WpbMLyVrx4WX59/f8MWMEwbI2Wqh9GVxWBFsVbCgc8NbAm5TDulv9omS1Dc62Gh0ZfkYiAiT2dnJzo7O+MeRlmmZXuBKRCinMBUKmQdTufw7oe+4zkTwzmrbCuBKkvFgclXwfKHp6jj/goXWwOJiIJsr0JUX3OnqLG4z4ybijOghRCwBcoGNsMXoLzQVqrdMzK8lWn3DFXrCi2ftveIhRFbwDQNBiKiyaQqMjpTMjpTiQm/l20LjBhWoCrlD0zpbHHVKp2zkM6YODCY9UKYG8jKLZ0rSU5roFa+KhVZtYoIWSpbA4mowbmBSOGy20SRJMldzr6+A18ul8Od33s27mFMCAMRNS1ZlrxKz0QJkf8/Jemcr0IVCkxuoApXsPb2jxRVuLJmZa2BgcCkqWhL+udd+QJX5Lyr/HGNrYFEFANvDpHCChERxYuBiKgKJEmCnsgvmzu9VZvw+4VbA4cywblUpULWoaEc9uSCbYRjaQ0MV6jCrYGjzbtiayARVcrmg1mJqE4wEBHVoWq3Bg4bY5t35Va43NZAfxvhWFsDo6pSkVWriJDF1kCiqctrmWOFiIhixkBENMXJsoQ2p7ozUW5r4FjnXQ1lgq2BbiCrpDVwtKpUYHn2MvOu2BpIVF8YiIioXjAQEVHF/K2BM6rwfoZlYzhi3tVoIevQUA67Q0uyp7NW2a/ltgaG5121jnHeFVsDiarDrpPnEBERMRARUWwSNWwNHG1xC681MGPiwECmqI1wLK2B5ZZcr2TxC7YGUjNihYiI6gUDERFNCdVuDcya9hifdxVcNdBf4RqtNTCZkCubd1XB4hZsDaRGYVn5qi4DERHFjYGIiChEkiQkEwqSVWwNLISpiJbAEiHr4FAWu0PHh3OjtwZWtFpgmXlXbghr0RSGK6oZVoiIqF4wEBER1VhCkTGtRca0luq2Bg6VC1mBila+NXD/kUxo+XYLVoWtgRUtblHmOFsDKcxmhYiI6gQDERFRA/G3BnZN8L3CrYHF8658c6284/l97/ePFI6PozWw7LwrTXHmWZWed8XWwMZncVEFIqoTDERERE2q1q2B4cCUzpkYyviqVk7IOjiUxe5Dw4EQNlprYEKRgqsFTmDeFVsD42FbFhSF33siih8DERERVUU9tAYOhVoD3eXaK2kNHG01wNGWZHfnXbE1sDKWbbNdjojqAgMRERHVnVq1BkaHq0Jr4JBvEYt01sR76ZzvIcNjaw0c67yrcOtgm65CU6duuLItC6rKjyFEFD/+S0RERFNaoDWwbeLvN1pr4FC4auW8/mAoi/85NBwIYWNqDZzgvKt6aw1khYiI6gUDERER0RhUvTUw51vEIpdfEdBbxCJnIp2xQqsDBlsDC4GsfGugLAEtWgVVK+e1mRnC/1jT8Mq+LGZk0tBVCboqI6FI0BXZ29YUCYo89qBlOXOIiIjixkBEREQUE1mW0JbMV3Oq1RoYWaEa5SHDH/paA90KVs60AZyADf/VD6C//J9DEkhIQEJ2fwnvtSoVtlXf8f5D/WhPzYD8+11IJmQkVQW683syIeereqoMPZHfTiUU6M4xnasMElEV1SwQ3XrrrfjP//xPbN++HZqmob+/v+icPXv2YNWqVdiwYQPa2trQ29uLtWvXsqeYiIhojPytgTOr0Bo4NDyCh3/57zj44RHYkgILEkwhwxRS/hdkWO5rIcOE5GwXzrHM/OuskJBG/nXhnBmQcjp2vLALGcMedW5WmK7KvuDkhinZ+x7kg5MSCFv5UCUHjqcSTuhy9keFMi7zTjS11Sx55HI5XHTRRVi6dCn+9V//tei4ZVk4//zz0d3djRdffBH79u3DV77yFSQSCXznO9+p1bCIiIioAm0tKXztq1+etK9n2wI5y8aIYSFr2MiYVj4oGRZGDBtZZztj+o7nfOe5xw0LGdPGcM7Ch2nDO3/E2Z91fs+NIYBJEgpByg1OvmDlBijdV9UqhK/C+SknrLn73fNTCSUQyhKKxABGNIlqFohuvvlmAMBDDz0Uefw3v/kNduzYgd/+9rfo6urC6aefjm9/+9u49tpr8U//9E/QNK1WQyMiIqI6I8sSknI+MEwG2863GI4Y+ZUDM77fi0KZsz//yx/KCvvTWQuH0jkvnI0YhfCVMSwYVun5XWGyhEDoCoexZKiqVbLlsFSVzB/iVAYwoth60zZu3IhTTz0VXV2Frunly5dj1apVeP3117Fo0aLI67LZLLLZrLc9MDBQ87ESERHR1CLLElKagpQ2OQHMskUgdHnVKl9VK3DcCJ3n329aGMqYOBjY716bD2RjCWCKLBWqWYFg5Yapwlwut0XRH7bcoOW/1tsXEcoSfFYX1ZnYAlFfX18gDAHwtvv6+kpet3btWq/6RERERNQIFFlyHuw7OV/PtOxQmHIqV7mIUGbayOSKQ5lbQcsYNgZGDBwYLN3CaJZZ4TBMlaWiKlUq1HJYWEgjokqm+c4rOW+scJwPS6bRjCkQXXfddbj99tvLnvPGG2/gxBNPnNCgyrn++uuxZs0ab3tgYABz586t2dcjIiIiajSqIqNNkdGmT87/+zYsO9RmmP/d337oP+5Wu/zByttnWOgfNtBnZJ32w+IWxnJLzIclFCmwQmF4zlbRghtFYWv0eWP+/fI4lqGneI3pb8nVV1+Nr371q2XPOf744yt6r+7ubmzevDmwb//+/d6xUnRdh65P0v9eISIiIqJRJZR8K1zbJDQfCSFgWCJioQ230hURtnIlQplTJTs0lAuGuFC1bAz5CwlFKp7zFTmXyw1f+f3e/C7fwhv+5eajQpmuygxgVTCm/2pnzZqFWbNmVeULL126FLfeeisOHDiA2bNnAwDWrVuHjo4OLFy4sCpfg4iIiIimFkmSoKkSNFVGe7L2X88NYNFzukILb5TY7y3g4ez/IJv1Vj+MamEUYwhgmiqXXGCj1EIaXosinwEGoIZziPbs2YPDhw9jz549sCwL27dvBwCccMIJaGtrw7nnnouFCxdi5cqVuOOOO9DX14dvfetbWL16NStARERERFQX/AGsYxK+nhACOdMe20IbEcHKP9fryIjhC2XF1bKxCD8DLKUKLK3R92Ky1CwQ3XjjjfjJT37ibburxm3YsAHnnHMOFEXBU089hVWrVmHp0qVobW1Fb28vbrnllloNiYiIiIiorkmSBN1pqetMJWr+9YQQgRUKC22G4YU2bG9OlzvXK2PayGZzwGs1H2ZNSUKMpShXfwYGBtDZ2YkDHxxER8dk5HYiIiIiIgKclkLDAAAkEonYW+oGBgYwe9ZMHDlypOJsENuy20RERERE1NgkSYKmaXEPY0K4MDsRERERETWthq8QuR1/g4ODMY+EiIiIiIji5GaCscwKavhA5P6hP3L8/JhHQkRERERE9WBwcBCdnZ0VndvwiyrYto29e/eivb099klcU8HAwADmzp2Ld999l4tUxIT3IH68B/HjPagPvA/x4z2IH+9B/MZyD4QQGBwcRE9PD2S5stlBDV8hkmUZc+bMiXsYU05HRwf/0seM9yB+vAfx4z2oD7wP8eM9iB/vQfwqvQeVVoZcXFSBiIiIiIiaFgMRERERERE1LQYiCtB1HTfddBN0XY97KE2L9yB+vAfx4z2oD7wP8eM9iB/vQfxqfQ8aflEFIiIiIiKi8WKFiIiIiIiImhYDERERERERNS0GIiIiIiIialoMRERERERE1LQYiIiIiIiIqGkxEJHn3nvvxXHHHYdkMoklS5Zg8+bNcQ9pylq7di0+/vGPo729HbNnz8bnPvc57Ny5M3BOJpPB6tWrMWPGDLS1teFv//ZvsX///phGPPXddtttkCQJV111lbeP92ByvP/++/jyl7+MGTNmIJVK4dRTT8XLL7/sHRdC4MYbb8TRRx+NVCqFZcuW4e23345xxFOLZVm44YYbMH/+fKRSKXzkIx/Bt7/9bfgXoeU9qK7nn38en/nMZ9DT0wNJkvCrX/0qcLyS7/fhw4exYsUKdHR0YNq0afi7v/s7DA0NTeKforGVuweGYeDaa6/FqaeeitbWVvT09OArX/kK9u7dG3gP3oOJGe3vgd/ll18OSZJw1113BfZX6x4wEBEA4Be/+AXWrFmDm266Cdu2bcNpp52G5cuX48CBA3EPbUp67rnnsHr1arz00ktYt24dDMPAueeei3Q67Z3zjW98A08++SQeffRRPPfcc9i7dy8uvPDCGEc9dW3ZsgU//OEP8dGPfjSwn/eg9j788EOcddZZSCQSePrpp7Fjxw5873vfw1FHHeWdc8cdd+Duu+/G/fffj02bNqG1tRXLly9HJpOJceRTx+2334777rsP//Iv/4I33ngDt99+O+644w7cc8893jm8B9WVTqdx2mmn4d577408Xsn3e8WKFXj99dexbt06PPXUU3j++edx2WWXTdYfoeGVuwfDw8PYtm0bbrjhBmzbtg2PPfYYdu7ciQsuuCBwHu/BxIz298D1+OOP46WXXkJPT0/RsardA0EkhDjzzDPF6tWrvW3LskRPT49Yu3ZtjKNqHgcOHBAAxHPPPSeEEKK/v18kEgnx6KOPeue88cYbAoDYuHFjXMOckgYHB8WCBQvEunXrxCc/+Ulx5ZVXCiF4DybLtddeK84+++ySx23bFt3d3eK73/2ut6+/v1/oui4efvjhyRjilHf++eeLSy+9NLDvwgsvFCtWrBBC8B7UGgDx+OOPe9uVfL937NghAIgtW7Z45zz99NNCkiTx/vvvT9rYp4rwPYiyefNmAUDs3r1bCMF7UG2l7sF7770njjnmGPHaa6+JY489Vnz/+9/3jlXzHrBCRMjlcti6dSuWLVvm7ZNlGcuWLcPGjRtjHFnzOHLkCABg+vTpAICtW7fCMIzAPTnxxBMxb9483pMqW716Nc4///zA9xrgPZgsTzzxBBYvXoyLLroIs2fPxqJFi/DAAw94x3ft2oW+vr7Afejs7MSSJUt4H6rkE5/4BNavX4+33noLAPCHP/wBL7zwAs477zwAvAeTrZLv98aNGzFt2jQsXrzYO2fZsmWQZRmbNm2a9DE3gyNHjkCSJEybNg0A78FksG0bK1euxDXXXIOTTz656Hg174E64dFSwzt48CAsy0JXV1dgf1dXF958882YRtU8bNvGVVddhbPOOgunnHIKAKCvrw+apnn/8Lq6urrQ19cXwyinpkceeQTbtm3Dli1bio7xHkyOd955B/fddx/WrFmDf/iHf8CWLVvw9a9/HZqmobe31/teR/37xPtQHddddx0GBgZw4oknQlEUWJaFW2+9FStWrAAA3oNJVsn3u6+vD7Nnzw4cV1UV06dP5z2pgUwmg2uvvRYXX3wxOjo6APAeTIbbb78dqqri61//euTxat4DBiKimK1evRqvvfYaXnjhhbiH0lTeffddXHnllVi3bh2SyWTcw2latm1j8eLF+M53vgMAWLRoEV577TXcf//96O3tjXl0zeGXv/wlfvazn+HnP/85Tj75ZGzfvh1XXXUVenp6eA+o6RmGgS984QsQQuC+++6LezhNY+vWrfjBD36Abdu2QZKkmn89tswRZs6cCUVRilbP2r9/P7q7u2MaVXO44oor8NRTT2HDhg2YM2eOt7+7uxu5XA79/f2B83lPqmfr1q04cOAAPvaxj0FVVaiqiueeew533303VFVFV1cX78EkOProo7Fw4cLAvpNOOgl79uwBAO97zX+faueaa67Bddddhy996Us49dRTsXLlSnzjG9/A2rVrAfAeTLZKvt/d3d1Fix6ZponDhw/znlSRG4Z2796NdevWedUhgPeg1n7/+9/jwIEDmDdvnvczevfu3bj66qtx3HHHAajuPWAgImiahjPOOAPr16/39tm2jfXr12Pp0qUxjmzqEkLgiiuuwOOPP45nn30W8+fPDxw/44wzkEgkAvdk586d2LNnD+9JlXz605/GH//4R2zfvt37tXjxYqxYscJ7zXtQe2eddVbRkvNvvfUWjj32WADA/Pnz0d3dHbgPAwMD2LRpE+9DlQwPD0OWgx8HFEWBbdsAeA8mWyXf76VLl6K/vx9bt271znn22Wdh2zaWLFky6WOeitww9Pbbb+O3v/0tZsyYETjOe1BbK1euxKuvvhr4Gd3T04NrrrkGzzzzDIAq34PxrQVBU80jjzwidF0XDz30kNixY4e47LLLxLRp00RfX1/cQ5uSVq1aJTo7O8Xvfvc7sW/fPu/X8PCwd87ll18u5s2bJ5599lnx8ssvi6VLl4qlS5fGOOqpz7/KnBC8B5Nh8+bNQlVVceutt4q3335b/OxnPxMtLS3ipz/9qXfObbfdJqZNmyb+4z/+Q7z66qvis5/9rJg/f74YGRmJceRTR29vrzjmmGPEU089JXbt2iUee+wxMXPmTPHNb37TO4f3oLoGBwfFK6+8Il555RUBQNx5553ilVde8VYwq+T7/Vd/9Vdi0aJFYtOmTeKFF14QCxYsEBdffHFcf6SGU+4e5HI5ccEFF4g5c+aI7du3B35OZ7NZ7z14DyZmtL8HYeFV5oSo3j1gICLPPffcI+bNmyc0TRNnnnmmeOmll+Ie0pQFIPLXj3/8Y++ckZER8fd///fiqKOOEi0tLeLzn/+82LdvX3yDbgLhQMR7MDmefPJJccoppwhd18WJJ54ofvSjHwWO27YtbrjhBtHV1SV0XRef/vSnxc6dO2Ma7dQzMDAgrrzySjFv3jyRTCbF8ccfL/7xH/8x8MGP96C6NmzYEPkzoLe3VwhR2ff70KFD4uKLLxZtbW2io6NDXHLJJWJwcDCGP01jKncPdu3aVfLn9IYNG7z34D2YmNH+HoRFBaJq3QNJCN+jqImIiIiIiJoI5xAREREREVHTYiAiIiIiIqKmxUBERERERERNi4GIiIiIiIiaFgMRERERERE1LQYiIiIiIiJqWgxERERERETUtBiIiIiIiIioaTEQERERERFR02IgIiIiIiKipsVARERERERETev/A27gky/qQSIAAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 1000x400 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lens = optic.Optic()\n",
    "\n",
    "# define beam expander properties\n",
    "radius_start = 5\n",
    "radius_end = 15\n",
    "\n",
    "# add surfaces\n",
    "lens.add_surface(index=0, thickness=np.inf)\n",
    "lens.add_surface(index=1, thickness=5, radius=-50, material=\"N-BK7\", is_stop=True)\n",
    "lens.add_surface(index=2, thickness=75, radius=np.inf)\n",
    "lens.add_surface(index=3, thickness=5, radius=np.inf, material=\"N-BK7\")\n",
    "lens.add_surface(index=4, thickness=50, radius=-300)\n",
    "lens.add_surface(index=5)\n",
    "\n",
    "# set aperture\n",
    "lens.set_aperture(aperture_type=\"EPD\", value=2 * radius_start)\n",
    "\n",
    "# add field\n",
    "lens.set_field_type(field_type=\"angle\")\n",
    "lens.add_field(y=0)\n",
    "\n",
    "# add wavelength\n",
    "lens.add_wavelength(value=0.633, is_primary=True)\n",
    "\n",
    "# draw lens\n",
    "lens.draw(num_rays=5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define optimization problem:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "problem = optimization.OptimizationProblem()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Add operands (targets for optimization):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Ray y-intersect at last lens surface\n",
    "input_data = {\n",
    "    \"optic\": lens,\n",
    "    \"surface_number\": 4,\n",
    "    \"Hx\": 0,\n",
    "    \"Hy\": 0,\n",
    "    \"Px\": 0,\n",
    "    \"Py\": 1,\n",
    "    \"wavelength\": 0.633,\n",
    "}\n",
    "problem.add_operand(\n",
    "    operand_type=\"real_y_intercept\",\n",
    "    target=radius_end,\n",
    "    weight=1,\n",
    "    input_data=input_data,\n",
    ")\n",
    "\n",
    "# Ray y-intersect at image surface\n",
    "input_data = {\n",
    "    \"optic\": lens,\n",
    "    \"surface_number\": 5,\n",
    "    \"Hx\": 0,\n",
    "    \"Hy\": 0,\n",
    "    \"Px\": 0,\n",
    "    \"Py\": 1,\n",
    "    \"wavelength\": 0.633,\n",
    "}\n",
    "problem.add_operand(\n",
    "    operand_type=\"real_y_intercept\",\n",
    "    target=radius_end,\n",
    "    weight=1,\n",
    "    input_data=input_data,\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define variables - let two radii of curvature vary:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "problem.add_variable(lens, \"radius\", surface_number=1)  # first surface of first lens\n",
    "problem.add_variable(lens, \"radius\", surface_number=4)  # second surface of second lens"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check initial merit function value and system properties:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "╒════╤════════════════════════╤═══════════════════╕\n",
      "│    │   Merit Function Value │   Improvement (%) │\n",
      "╞════╪════════════════════════╪═══════════════════╡\n",
      "│  0 │                 49.144 │                 0 │\n",
      "╘════╧════════════════════════╧═══════════════════╛\n",
      "╒════╤══════════════════╤══════════╤══════════╤══════════╤══════════╤════════════════════╕\n",
      "│    │ Operand Type     │   Target │   Weight │    Value │    Delta │   Contribution (%) │\n",
      "╞════╪══════════════════╪══════════╪══════════╪══════════╪══════════╪════════════════════╡\n",
      "│  0 │ real y intercept │       15 │        1 │  9.22613 │ -5.77387 │            67.8366 │\n",
      "│  1 │ real y intercept │       15 │        1 │ 11.0243  │ -3.97572 │            32.1634 │\n",
      "╘════╧══════════════════╧══════════╧══════════╧══════════╧══════════╧════════════════════╛\n",
      "╒════╤═════════════════╤═══════════╤═════════╤══════════════╤══════════════╕\n",
      "│    │ Variable Type   │   Surface │   Value │ Min. Bound   │ Max. Bound   │\n",
      "╞════╪═════════════════╪═══════════╪═════════╪══════════════╪══════════════╡\n",
      "│  0 │ radius          │         1 │     -50 │              │              │\n",
      "│  1 │ radius          │         4 │    -300 │              │              │\n",
      "╘════╧═════════════════╧═══════════╧═════════╧══════════════╧══════════════╛\n"
     ]
    }
   ],
   "source": [
    "problem.info()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define optimizer:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "optimizer = optimization.OptimizerGeneric(problem)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Run optimization:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "  message: CONVERGENCE: REL_REDUCTION_OF_F_<=_FACTR*EPSMCH\n",
       "  success: True\n",
       "   status: 0\n",
       "      fun: 2.2236952958197887e-12\n",
       "        x: [-1.214e+00 -1.641e+00]\n",
       "      nit: 19\n",
       "      jac: [-1.639e-04  2.945e-05]\n",
       "     nfev: 111\n",
       "     njev: 37\n",
       " hess_inv: <2x2 LbfgsInvHessProduct with dtype=float64>"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "optimizer.optimize(tol=1e-6)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Print merit function value and system properties after optimization:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "╒════╤════════════════════════╤═══════════════════╕\n",
      "│    │   Merit Function Value │   Improvement (%) │\n",
      "╞════╪════════════════════════╪═══════════════════╡\n",
      "│  0 │            2.51819e-12 │               100 │\n",
      "╘════╧════════════════════════╧═══════════════════╛\n",
      "╒════╤══════════════════╤══════════╤══════════╤═════════╤══════════════╤════════════════════╕\n",
      "│    │ Operand Type     │   Target │   Weight │   Value │        Delta │   Contribution (%) │\n",
      "╞════╪══════════════════╪══════════╪══════════╪═════════╪══════════════╪════════════════════╡\n",
      "│  0 │ real y intercept │       15 │        1 │      15 │ -7.01452e-07 │            19.5393 │\n",
      "│  1 │ real y intercept │       15 │        1 │      15 │ -1.42343e-06 │            80.4607 │\n",
      "╘════╧══════════════════╧══════════╧══════════╧═════════╧══════════════╧════════════════════╛\n",
      "╒════╤═════════════════╤═══════════╤══════════╤══════════════╤══════════════╕\n",
      "│    │ Variable Type   │   Surface │    Value │ Min. Bound   │ Max. Bound   │\n",
      "╞════╪═════════════════╪═══════════╪══════════╪══════════════╪══════════════╡\n",
      "│  0 │ radius          │         1 │ -21.3565 │              │              │\n",
      "│  1 │ radius          │         4 │ -64.0695 │              │              │\n",
      "╘════╧═════════════════╧═══════════╧══════════╧══════════════╧══════════════╛\n"
     ]
    }
   ],
   "source": [
    "problem.info()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Draw final lens:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA0QAAADMCAYAAABTEmMGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/TGe4hAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA/LklEQVR4nO3de5Qbd303/vfcNDO67329u7bjGJNgJ3FCQoIJT4EHl5TSAk0PLTxpcIFTHtLkQDA/CGkJKdBgEk4phXJCy1MaekoLpQ/hkhJKcEIgT5w4sWNyNyYxdry72rsuK2lGl/n+/pA0klZar7230a7er3P27KxmJH21Y+/qvZ/5fr6SEEKAiIiIiIioDcleD4CIiIiIiMgrDERERERERNS2GIiIiIiIiKhtMRAREREREVHbYiAiIiIiIqK2xUBERERERERti4GIiIiIiIjaFgMRERERERG1LdXrASyV4zgYGRlBKBSCJEleD4eIiIiIiDwihEAqlcLAwABk+cxqP2s+EI2MjGDjxo1eD4OIiIiIiFrESy+9hKGhoTM6ds0HolAoBAB44cXj7jYREREREbWfVCqFreduOatcsOYDUeUyuVAohHA47PFoiIiIiIjIa2czlYZNFYiIiIiIqG0xEBERERERUdta85fMERERERGRN4QQyOfzAABN09Zk12cGIiIiImrq2WeeQSKRWPC4giOQypU+ZvMCyfI2VB19fX2rMFIi8oooFhB/5D8BAHs/8v/B5/N5PKKzx0BEREREdYQQOPbSKL7+vftQUEzY0JAVCixHRVaosISCrFP+LFTkhNLwGAocaFIG/hctD14BEa0WRRTxe14PYokYiIiIiNqMEAKTszmMxLMYjls4Fc+628Pl7WzeAbDdvY9fkxExFUQNFf2GgoihIGqq7ueooSBS+dpQEJ+IYf/+/Xjfe/8MPT093r1YIlpRuVwOX/ibg14PY0kYiIiIiNYZxxGYmLXrAs5w3MLwTBbDiSxG4hbsguMeHzZUDEZNDEQNXLm1C4NRE35YOPjAvXjLG38L52zohk89uz5MWm8vJEnC8KlTDERE1NIYiIiIiNaYoiMwnrIxHM+WPmYsjCSypcATL23ni8I9PmpqGIwaGOww8bptPRjsMDAYNcsfBkKG1vAcsVgMLz6YQZdfOeswBJQmV3d0duLUqZdw8SWXLOn1EhGtJAYiIiKiFlMoOhhL2eWAU72UrVTtsTCasFBwqoGnM6CVKzwm3rghhKHy9mDUwEDURFD35td9T08PTp065clzExGdKQYiIiKiVZYvOoglrIagU9mOJW0UawJPd9DnVnQuGoy4YWeow8SGiAG/rzV/nff09ODo888jk0nD7w94PRwioqZa8ycoERHRGpYrOBhNVAOOG3xmshhJWBhLWqjJO+gN6e7la5dsjGKww8RAxHA/G1pjF7e1oDJ3aHh4BNu2bfN4NEREzTEQERERnSUrX8RIwmpoVjA8Uwo+E7M2RDnwSBLQVwk8HSYu39JZms9TrvhsiBiLmqOzFgQCAZimiZGRYQYiImpZDERERERzZHNFjCSyOFVpUjCnJfXEbM49VpEl9IdLgWdzl4nXbO10O7YNRU30hddv4FmIJEno6urCyPCI10MhIpoXAxEREbWdWbuAkZp5O6fmzOGZTufdY1VZwoaIgYGoga09AfzWtu5qhafDRF9Ih6q0Z+A5E93d3XjuuecghIAkSV4Ph4ioAQMRERGtOykr37xhwUwWwwkL8Uw18GiKhIFIaf7OeX0h/M/zekttqSOlwNMb0qHIfCO/WF3d3bBtGzPT0+js6vJ6OEREDRiIiIhoTRFCIGkV3KBzym1NXQ0+SavgHu9T5VJFJ2LigsEIrtrRV+7SVgpBPUEdMgPPiunq7AQAjMZiDERE1JIYiIiIqKUIITCTyWMkYbnr8IzELZyqCT2zdjXwGJrsNijYuTGC372wH0M1a/B0BXwMPB7SDQPBYBCx2Ch27Njh9XCIiBowEBER0aoSQmA6nau7pK3SqW0kUdrO5Iru8X6fUl5o1MCrNndgcGdpe7DDxFDURIdf49yUFtfZ2YlYLOb1MIiImmIgIiKiZeU4ApPpXLkzW7kVdaI8fyduYSSRhZV33OODulq6pK3DxKu3dLnblcVHoyYDz1rX2dnJxgpE1LIYiIiI6Kw4jsDErD2nJXVNW+qEhVyhGngiplpuQ23if2zrci9vq3RqC5uah6+GVkNnZyds20YikUA0GvV6OEREdRiIiIioTtERGEtaDfN2Kpe3jSYs5IvCPT7q1zAYNTEUNfCG83qqYafDxGDERNDgr5p211FurDA2NsZAREQth7+liIjaTKHoIJa0GzqznZopVXdiCQsFpxp4ugI+t0HB9g3hupbUAxEDAZ2/Suj0TNOEYRgYHx/Deeed5/VwiIjq8LcYEdE6kys4iCWtxmYFM6XbYkkLNXkHPUGfu8jozqFIqbJTvqxtIGLA9CnevRhaFyRJQjTagYnxca+HQkTUgIGIiGiNyRWcuoBTN38nnsVYyoYoBx5JAnpDenkOj4FLN0fd7Urg0TUGHlp50Y4oYqOjXg+DiKjBigain//85/j85z+PQ4cOYXR0FHfffTfe/va3u/uFELj11lvxta99DfF4HFdeeSXuvPNObNu2bSWHRUTU0qx8saFZQWX7VNzCRMp2j5UloC9sYDBqYGOHiSu2dGIwapTW4ekw0R824FNlD18NUUk0GsXzzz2HfD4PTWMjDSJqHSsaiNLpNHbu3In3vve9uPrqqxv233HHHfjSl76Eb3zjG9iyZQtuueUWXHXVVXj22WdhGMZKDo2IyDOZXKE6b6e2YUF5Ds/kbM49VpEl9IdLFZ5zugK4cmtXee6OicEOA/1hA5rCwEOtr9JMYWpqCv39/d4OhoioxooGoje/+c1485vf3HSfEAJf/OIX8YlPfAJve9vbAAD/8i//gr6+Pnzve9/DO9/5zpUcGhHRipm1C/Wd2WayGE5Y5XV4spjJ5N1jNUVCf7h0+dq2vhBeX9ulLWqiN6RDZeChdSASiQAAJicnGYiIqKV4Nofo+PHjiMVi2L17t3tbJBLBFVdcgQMHDjAQEVHLSmbzGElY5XV46ltSj8QtxLP1gafSoOAVG0LY/Ypedw7PUNRET0iHInOhSlr/fD4fAoEApiYnvR4KEVEdzwJRLBYDAPT19dXd3tfX5+5rxrZt2Hb1+vlkMrkyAySitiSEQCJbcBsUDDdZiydlFdzjdVV2GxRcOBjBm3f0u18PdpjoDvggM/AQAQDC4TAmGYiIqMWsuS5z+/btw6c+9Smvh0FEa5QQAjOZfLU724xV7thW7tSWyCJtF93jTU12KzyXbIzg9y7sLy84Wlp4tCvggyQx8BCdiXAkwtbbRNRyPAtEleuHx8bGsGHDBvf2sbExXHzxxfPe7+abb8bevXvdr5PJJDZu3Lhi4ySitUUIgal0rmHeTm2ntmzecY/3+xQMlQPO5Vs6MBgdqFZ4oiY6/BoDD9EyiYTDOParX8FxHMgy58YRUWvwLBBt2bIF/f392L9/vxuAkskkHn30UVx33XXz3k/Xdei6vkqjJKJW4zgCE7N23byd4TlzeOxCNfCEDNVtUrBra1epHXXUwEA58ERMlYGHaJWEw2E4joN4PI7Ozk6vh0NEBGCFA9Hs7Cx+/etfu18fP34cR44cQWdnJzZt2oQbb7wRf/3Xf41t27a5bbcHBgbq1ioiovZSdAQmUnbDvJ1KS+rheBb5onCPj5qaW9F53bYeDHYYGIiYGOooLToaNrneCVGrCIXDAICZ6WkGIiJqGSsaiB5//HG84Q1vcL+uXOq2Z88e3HXXXfjYxz6GdDqN97///YjH43jta1+LH//4x1yDiGgdKxQdjKXsasOCStApX9oWS1p1gafDr7kNCv5nf6g0f6fcpW0wYiJorLmpkERtKxAIQFEUzMzMeD0UIiLXir6TeP3rXw8hxLz7JUnCpz/9aXz6059eyWEQ0SrKFx3EkhaGZ6y6Tm2VS9piSRtFp/pzoTvocwPOBYPhujV4BqIG/D4GHqL1QpIkBIMhBiIiail8p0FEZyVXcDCaKDUoODWnWcFw3MJY0kJN3kFPSMdQOeBcvDFaE3ZKl7SZPsW7F0NEqy4YCmImzkBERK2DgYiI6tj5ojtXpxJ2TtW0px5P2agUfiUJ6AvpbsC5bHNHee6OicEOAxvCBnSNgYeIqkKhEMbHxrweBhGRi4GIqM1kc8XSujtzmxWU5/BMzObcY2UJ6A+XKjqbu0zsOrezPIfHwFDURF/YgE9l61wiOnPBYBDHfvUrCCHY4ZGIWgIDEdE6k7YLdQFnbkvqqXQ18KiyhP6IgcGogXN7Avgf27rrWlL3hXVoCgMPES2fYDCIYrGIdDqNYDDo9XCIiBiIiNaaWauAU3OaFZQuayttxzN591hNkbAhUqrwnNcXwhvO68Fg1CytxdNhojekQ5H5F1oiWj2BQAAAkEgkGIiIqCUwEBG1ECEEklZhTkvqcvApbyeyBfd4nypjIFK6fG3HhhDe9IpeDJbX3xnqMNET1CEz8BDRMik6AlbBgV0ofbbyDuyCA6sgyp/rt+2CgJV36u6TyRVw0n45Hr/71/AHRpo+z2ka1J6Rhe4vcPoDFr7/0p5/oUc43f0Xfu4FHnvB+y+wf4H7L/QAK/38S339Cx2wlH87S3/tS/t3u5DFPr8sinhD0z1rBwMR0SoSQiCezVc7s82Ug09NxWfWrgYeQ5NLDQqiBnZujOB3L+x3u7QNRk10BXwMPES0oKIjMGsXkLTymLUKSFoFpKwCUlYeKbuAZLb8tV2+zSpgOpXFqHUBvnfPBPLOBKyCg4JzZs+nyIChyjBUGboq1W3rqgxTdhBQi+gO6vM/yAI/2qQFDlhoetJCPzkXvv/Snn8hp5tftfSxL/65z8RKP/9Kn7sl7j7t86/0v1sv/t+IYgHFI6e/X6tjICJaRkIITKdzGK6ZwzNSM4dnOG4hkyu6x/t9ihtwLt3cgbfurIadwaiBzoCPk46J2pwQpYpLXYixqgEmmc2Xw0594Jktf52yC0jbxXkfX1dlhAwVIUNF2NAQ0lVE/Rq6DECdOY7tW7rQEQ6UA40MoxxqDK1mW63fVpXT/9z6wfefwMvOjeK333TRcn+7iGiV5XI5fOGI16NYGgYiorMghMDkbM4NOqcqQafSqS2ehZWv/gk1oCtuwHn1lq5qw4Jyp7aoqTHwEK1zlepMbZhJWgW3YuOGmzn7S4Gn9HW+2PxSFUkCgrqKsKGWP2sImSo2dpgIGxqCRmlfJfCE9NL+kF65XZu3U2QsFsNdx3+C392xE11dXcv6PfH7A0ilksv6mEREi8VARFTDcQQmZu2mLalPxbMYTViwa64ZCRtqKfB0mLhya1epYUGHiYFy1SdsqAw8RGucnS+Wqy8FzFr5UlixG6s1zfcX6i6DncunyjVhphRQwqaGoQ6zHGa0mjCjNoScgE9dk5fN+v0mUqmU18MgIgLAQERtpugIjKdsN+CMzGlJPZLI1v0lNurXMBgxMNhhuh3aBtw5PAZChubhqyGihTiV6kxNgEnWVGTq5tTYBaSy+bowk7TyC1ZnQnpNBcbQMNhhuuEmpKsImzXVm7pqjdq2Cxf7/X6McXFWImoRDES0rhSKDmJJu9ygoH4NnuG4hVjCQsGpvrnpDGjuJW2v2BAqtaMuh56BqImgzv8iRF6y88XSPJnaOTFzAk1qTpipHltAOleYt3OSpkh1IaV0uVkp0JRCjlYTZqqBJ1Su0KzV6kwrMP1+pNNpLs5KRC2B7/ZoTckVHIwlK3N3rGrwmSlVd2JJG8WawNMT9JUDjomLBiN1Lak3RAz4ffwvQLRSHEcgnSs0NANw589kq3NkGgJOOeTkTtPWzJ07U7mETFcxEDVxXk3FpvYSs7rLzdq4OtMKTNOEEAKZTMZdl4iIyCt8N0gtJVdwMJKoXsp2qmYOz3A8i7GU7f61V5KAnqBensNj4JWbom6zgoFoKfgYfMNDtGi5gtO0IlMNM3NCzpyKzax95tWZkK4iZGrYEDFqQs6c+TOm5m4HdJWLCq9hpmkCANLpNAMREXmOgYjgOA4ymcyKP48QApmcg9mijJGE3bQl9XjKdo+XJaAvbGCgPIfn8i2ddWvwbIgY83ZHImp3jiOQyRXdtsupbLUFc7P5M7Utmiv77QWqM/WT/VX0Rwy8vC9YN5+mdn/Q0MqVGlZn2p1hGABKgYiIyGsMRIR77/0RnnryybO+nxCADQWW0GAJFRbU0mehwZ7ztSVU2FBRRDXAKLKE/nCpwnNOVwCvKXdpqzQs6Asz8FD7yhWc6qKZTebPNISYORWb1ALVGfdyspr5Mf0Ro27S/9z9lYATZHWGloiBiIhaCQMRYWJ8AgMDA9j68vOQyom6j6QtkMo5jbfnBGZzAnPfb0kAQj7J/egqfw7rEkI+GWMnfo2ekI73vetq9IV0qAoDD60/Qgik7UozgPnXmXHXoSmHmdqKzemqMwFdcVswh2rCzDajZj0aU2voflZpBqCrMieyk6c0TYOiKKtydQIR0UIYiJbRXQ+fwF0HTng9jLM2O9uHvFBgvdC4JoSmSIgaCiKmioihYDCoIGKoiJoKooaKsKG42xFTQUhXIJ/mjdYBOYZ0JoPBqLmSL4loSXIFp+mimQ0tmudUbCpfz9oFOPNUZ1RZapjsHzJU9FXmzpTbNFe6nFVCTLCmYsPqDK0HhmEgy0BERC2AgWgZvWJDCFdfMuD1MM7a0aNHkc/M4IpLLiiFn3IAihoKTI1/Saa1RQiBdK7oXk5W24I5ZedrLi2bP+RY+fmrM36fUn8Jma6iN6xja0+goYVzbeCpHGvw/xQRgFIgymQZiIjIewxEy+iKLZ24Ykun18M4a4cicez/6VPYvfXVUFX+kyBv5YtOXdvl2rVlGteZqW/RXKnQLFSdqQ0oYVNFb0ifE2ZqGgK4jQM0BHWFl3kSLROfT0c2m/V6GEREDEQEDA4MwnEcTE9Po7e31+vh0BpW6iRYdANKMltu0VwzfyZpNc6nqYacPLILVGfcS8j0Usey3lCpOhNqFmjKgSdYnm/D6gxR69B1HywGIiJqAQxEhJ7eXiiKgqmpKQaiNlepzlRDTGNDgGRNRaY28MyWbyvOU55RZaluUcxKgOkJ6fUdzWrmz9RuszpDtL74fD7MxONeD4OIiIGIAEVR0NPTi+mpKa+HQktQqc7Mrcg0bdHctGJzZtWZ2jDTG9JxbmXujK4i5AaYcqCpueTM1BRWZ4jI5dN12Jbl9TCIiBiIqKSvvw8nT6y9DnnrSaHoIGVXGwDUTfavhJvTNgOYvzqjyJIbVGo7lnV3BxpCTu3+SpgJ6io0VmeIaBn5NA22bS98IBHRCmMgIgBAT08PnnrySTiOA1nmG9+zJYRANl9sbNHsNgNoXGdmbsjJ5IrzPr6pyXWBJWRo6A7q2NIdcNejaehoVrMmjd/H6gwRtRbN54Nt2xBC8OcTEXmKgYgAlAKR4zhIpVKIRCJeD2fVFYoOZu2iG1CS2ULDOjRNF9ismT9TmKc6I0uoCyyVRTO3dPvLk/2bdDQrX24WLFd1WJ0hovVG0zQAQC6Xg67rHo+GiNoZAxEBALq6ugEAiURizQUiIQSsvFOtyLjrzMwzfyZbvs2uzp9ZqDoTrMyHKQeUzoAPm7v8DS2aawNP5VhWZ4iIGjEQEVGrYCAiAEAgEIDP50MqmVz15y46oqEa07RF83whZ4HqTO2aM5UKzDld/vqOZ6ZWt117LKszRETLrzYQERF5iYGIAACSJCESjSI1O9t0vxACdkHALjiw3M+lD7sg6j5beafuOLsgkC2UbotNBJAtmPjPv3/Yrdik7fmrM0Zl7oxenezf4fdhU6e/2gzALC2aWdvRrHKJWoDVGSKillRZCDyfz3s8EiJqdwxEy+h7R0Zw9xMjXg9j0canBpAfK+Ku4RdhFxxk3UBTCjdnylAlGKoMXZVhqFLps1ba9skChlrExed21oWc8Jw5NJU1aXwqqzNEROsRAxERtQoGomVkagq6gj6vh7Fo+VkV6VQG23pM6OUgY2qyu10KOaWwU7etVbd9inTaisyBh4eRzmSw53fPX8VXRkRErUZVFAAMRETkPQaiZXTVjj5ctaPP62Es2qOPZPDQQ7/EO3/r1V4PhYiI1jmlXCEqFBiIiMhbvB6JXIFgEPl8HoVCweuhEBHROqeUK0SFPH/nEJG3GIjIFQgEAADZbNbjkRAR0XrnBqLi/I11iIhWAwMRufx+PwDAtm2PR0JEROtdZb5pscgKERF5i4GIXIZhAOCaEEREtPIkSYKiKCiyQkREHmMgIldlpfA8AxEREa0CWZbhzLOwNhHRamEgIpfPV2oZzhaoRES0GiRJgnAcr4dBRG2uJQLRV77yFZxzzjkwDANXXHEFDh486PWQ2pIsy1AUhRNciYhoVUiSBIeBiIg85nkg+va3v429e/fi1ltvxeHDh7Fz505cddVVGB8f93pobUlRFDgMREREtAokSYIjGIiIyFueB6IvfOEL+LM/+zO85z3vwfbt2/HVr34Vfr8fX//6170eWluSZRlF/rWOiIhWgSRJAKcQEZHHVC+fPJfL4dChQ7j55pvd22RZxu7du3HgwIGm97Ftu64tdDKZXPFxnqlsrogXJ9NeD2NJphw/TqYE9ImVWYtoOKvAsjU8M9I6542IiBpNTWUw6fjx4kweU87K/E6YKBj4TbKIKH8nEK1ZhZq551a+iPKU9DXF00A0OTmJYrGIvr6+utv7+vrw/PPPN73Pvn378KlPfWo1hnfWXpxM4+qvPuL1MJboXOApAE8dX6HHDwEI4etr/vtERNQOtuOH+6cBTK/Q428FHrGBR/g7gWitUlHEtWZp+/hkBjsDprcDWgRPA9Fi3Hzzzdi7d6/7dTKZxMaNGz0cUdW53QF89wOv9noYS/Lv//ZNbD7nHLzsZS9bkcd/8sknYVkWfu/3fn9FHp+IiJbH1NQkfvCDH+C1r30tIpHIijzH/p/+FNu378DOiy9ekccnopVXyOdx7789AQDY0u33eDSL42kg6u7uhqIoGBsbq7t9bGwM/f39Te+j67q7Xk6rMX0KdgyEvR7GknTJGWwOSXh5z8qk+ymziLTIr/nvExHReheTM+iWMzi3Q0NX18r8TvilksXmiMzfCURrWC6Xw73lbUNTPB3LYnnaVMHn8+HSSy/F/v373dscx8H+/fuxa9cuD0fWvgqFAmRlbf5jJiKitUeWPO/vRERtzvNL5vbu3Ys9e/bgsssuw+WXX44vfvGLSKfTeM973uP10NqOEALFYhEqAxEREa0CRwhIsuT1MIiozXkeiP74j/8YExMT+OQnP4lYLIaLL74YP/7xjxsaLdDKy+VyAABV0zweCRERtQPhOJBlVoiIyFueByIAuOGGG3DDDTd4PYy2V2ln7mMgIiKiVeA4Di+ZIyLP8acQuWzbAgD4WrRpBRERrS+O40DhZdpE5DEGInJlMqWF91q1ix8REa0fjuNACMFGPkTkOQYicmXSaQCAYRgej4SIiNa7YrEIAFDVlrh6n4jaGAMRuWbTs1AUBRrnEBER0QpzGIiIqEXwp9AyeuJkHL88lfB6GIt27FgSExjAvUfj0FUZhipDVyUYc7e10rYmS5AktkslIqKzV2AgIqIWwZ9Cy+jgb6bxD7847vUwFq2QL6IoevHgg6NndLwsoRycpHJgagxNldsrx4yN60ARiP5yBCFdRdjQEDRUhA0VIUNFwKdC5poURETrXuWSOU3jWxEi8hZ/Ci2j//1b5+J//9a5Xg9j0b5x1z/D7/fj1a95DeyCgF1wYBUc2AUBq+DAytds19xul/dZ5fvYBYFswUEm52AmU4BVc9xsVodVNPCT//t00zFIEhDUKwFJQ0hXETLVpuEpZKgI6VrN/tJ9fCqvBCUianWFQgEAoGk+j0dCRO2OgYhc8XgcvX19kCUJpibB1JY/WBx4+GGkMxm8839di6RVQMouIJXNI2UXSl9bBaSsfPlz6bZZK49TM1mk7JR7+6xdmPc5dFVG2FDL4ak+VIUMDaFyqJq7v7Lt9ymsUhERrbBCPg8AnLdKRJ5jICIAgJXNIpvNIhQKrcrz6ZqCHk1BT2hxLb6LjkC6HKJm7QKS2Xx1uyZQVQJWIpPHS9PZuv35omj62HK5ShUqV5yq4alckaqpUIWNxq+DusoqFRHRAvLlCpHPxwoREXmLgYgAAFPT0wCASDjs8UjOjCJLCJsawubi/rIoROnyvkoVqrYqlSoHrNk5VauXprNIWalyNSuPtF2c9/ENTa6rSM0NTGGjUrWaE7jKxwd8ChtWENG6VqkQMRARkdcYiAgAMDk5AQAIRyIej2R1SJIEQ1NgLLFKVak4zVpzLvmzC0hmC26omrXyiGcqoaqAlL1wlSpkqOXwVB+qKlWr0qV+mjt/KlhTsWKViohaXY6BiIhaBAMRAQAmxicQDofZ/vQsKLKEiKkhsoQqlZV3SvOoyhWqUniac9lfTeg6OZ1Fykq6x2ZyC1SpjOaBqek8qtoqFqtURLTC8vk8NE2DLPOPN0TkLb77JQDA2FgM0Y4Or4fRViRJgulTYPoU9C6ySlUoOpi1i27FqXbeVLL263Lomk7ncGI6U3eJYMFpXqVSZAlBXal2+zNUhM3m4WluxaqyX1P4RoeImsvnctD1xf3sIyJaTgxEBCEExsbGsH3HDq+HQmdJVWRE/TKi/qVVqZJWec5UtrBgxeo3sxnM2tXAdboqlVmpUlUqUnopVFVbq6tN91dClZ9VKqJ1y87lYBiG18MgImIgImB6ehq5XA5dXV1eD4VWWW2Vqm+Rj+FWqRrmTuXrtiuBamo2h99Mpt1W67MLVKlqK051nf/0+kYUDXOsylUrlVUqopaUYyAiohbBQEQYGR4GAHR3d3s8ElqLlqNKlc0X6y75cwOUXUAqWz+PKmkVcHwy4wawhapUfp/StCJVCUyNC/7WBCpDhamxSkW0EnK2Db/f7/UwiIgYiJbT8ck0TkxlvB7GWXv86WGkzA0YnnUQMQqIGAoULkxKq0SSJPh9Kvw+FX2L7PqeLzqYtQtuYJqv818lcE3O2vjNZLq6OLBVQPF0VSo3PDU2p5i74O/cxX9ZpSJqzrZtdPLKBCJqAQxEy+i/norhyw+84PUwFkEBMIj//M6LAAAJQEhXEDEVRA0FYUNF1FAQNVVEaj5XtqOGwhbP5ClNkdHh96HDv7j2vUIIZHLFhrWnGuZR1YaqyXTdOlbZvDPv4/t9StPwVNflz6yZRzXnMkBWqWg9sm0bfpMVIiLyHgPRMnrPazbjjy8b8noYZ+3b3/oWFN3E0LbtSFhFJO0iEpaDhFUsfdgFHE3ZSFoOEnYRdqHxL+mmKiFsKIjoCiKG7Aam8JyvR5J5mHxfRy1GkiQEdBUBfWlVqpRVqAlVzTv/VfZPzNp4sVylqgSueYpUUGWpvptfkwV/Gxb/NavHBnWFVSpqOZZl8ZI5ImoJDETLqPKGaq3pCek4fvzXSAz/2r1NAhAtf9TRgIIqwxIqLKilz0KFJTRYlgorqyImVJyA5u7L1f0z6wEA/NPnHsBg1MRg1Ch97jAxEDEw2GFiMGoiuAa/j9TeNEVGZ8CHzsDSqlTu2lPZ/PwVq/L+iZRd05zizKpUpSpUeZ6UOf88qlA5VFUqVoYms0pFy6ZYLCKfzzMQEVFL4LtOwu+/9a0YHj61Yo9fcASSVhFJu1R1EkYEExkHw/EsRuJZ7H9+AiOJLPLF6p/Ho6aGgUpYipru9lA5OIUXuRgqUauqrVL1L/IxKlWqakOKfF03v9puf0krj/GkjResdF2r9dNVqZp28mtY2LdJ5z+99DXnJlKFZVkAwEBERC2BgYjg9/uxbdvLPR2D4whMzNoYjltuUBqOWxieyeLBYxMYjlvIFap//Q4Zan1QKlebBsoBKmKq/Gs2tZ3lqFKlc8WG8DS3y9/cUFVtTpGHtUCVKjxnQd+6eVTNFvytCV2sUq0fVjYLAAgEAh6PhIiIgYhahCxL6Asb6AsbeOWmaMN+xxGYSufcoHQqni0HJwsPvzCF4Xi27o1YQFcwGClfktdhYiBiYrCjWnHq8Gt8Y0U0hyRJCOqly+T6I4tbHyZXKHX8qwQmN0SVK1a14SllFUpVqvG0uzjw6apUmiLVVZyqzSlq1qYy529eEdRZpWoV2XKFKBAMejwSIiIGIlojZFlCT0hHT0jHzo2N+4UQmMnkS0FpphSUhsuh6eDxGZyKj9StVWNqcrnCVDOPKVoNUF0BHwMT0SL4VBmd6tKrVLWNKFLzzaMqfx1L2nXt1U9XpQroSkNFqi486SrCZuP+SrVKV1mlWg7ZcoWIl8wRUStgIKJ1QZIk91KhiwYjDfuFEIhn89WgNJPFcKK0ffhkHD98MoZZu+Aer6ty3RymwTnzmXqCOmT+pZlo2dVWqTYssUqVrHT2m3sJoHspYGn/WNLCsfGazoB2AWKBKpXbza/Zgr+1nf/KAatybIBVKgClQOT3+6EoitdDISJiIKL2IEmSu07NjoHmfZWT2Xy5qlStLg3HLTw1nMCPn4khka0GJk2RSuGopjNe7XymnpDONz1EHllqlcpxSh3/3LWnsvPNo6rujyXtuv12Yf4qVbC2+qTXX/bXdB6VriJUrlqFDRW6tvZDRCaTQZCXyxFRi2AgIioLmxrCpoZXbGgemGatAoYT5epS3CrPZ8ri2dEUfvrcOGYyefdYTZHQH65vKT7UUW380BfSuS4MUYuSy+s+BY3F/4rMFZzqJX9uqKrt8le/4G8sUa1SVS4JPF2Vqm7x3jmBqbFaVXNbufrmdYU7k8kgFFrkol9ERMuMgYjoDAUNFecZIZzXF2q6P20XMJKwypfjZTE8U6o0HRtL4WdHJzCVzrnHKrKE/rDeMHep0vyhP2xAY2AiWrN8qoyuoI6uoL6o+zuOQDrXPDzVdvmrrViNxrN1i//m5qlSSRIQ8NU2nKiGKncelaHByaVxvNCBQ6M2eqw0DE2GoUrQVRmGWt1eTDU8m8mgp6dnUd8bIqLlxkBEtEwCuoptvUFs621+GUg2VywFppoOecPxLI5PpvHQrycxMVsNTLIE9IUNd+5SpfnDUHl7Q8SAT2VgIlqvZFkqV3U0DCzyMex8sS5AuSHKLiCZrXb7qwSqkXgWR2tarc9aBQhsxc/+XxxAfN7n0WQJuiqVA5MMvRyUzJptQ5VhaKWvDVXG81MmhoMmYg+fWOSrKwW7Je3H6Q9Y+P4LWMHnX+i5F2r8sfD9FzjA4+df6utfyEr+21n6a1/aN2clnr9YyDc5cm1hICJaJaZPwdaeALb2NF93w85XA1MlLJ2KZ/HSTBaPHJ/GeMp2L6GRJKA3pLtzl4bmLGA7EDHWxTwDIlo8XVOgawq6F1mlGhkdxf/553/B69742zBDUVgFB3bBgVUQsPLV7dJnp7xflLbz1e14tlB/XN5ByurBkRfykI//elFjE/NdT1jZv+D9F9i/xOdfyFKef6VfO9HZUlHEtabXo1gaBiKiFqFrCrZ0B7Clu3lgyhUcxJKlS/JOxetbix8+EUcsadWt39IT9NXNYao0fxiMGhiImDB9DExEND9ZkuCTiugwFXRFFtegoplEIoEffP/7+F/XXINNmzYv2+PS6lgwkK1w2FzpwLeSzy8WuPdSv3cLHbBSz5/P5/B/vvLE6e/c4hiIiNYInypjU6cfmzqbr9uRL1YCk4WRmuYPw/EsjrwURyxpo1iTmLoCPrfJQyk4GaXFbMsBKqDzxwMRLb/Z2VkAQDjcuEQCtb4FL4dbcr8Odmhda3LK2i878h0P0TqhKTI2dvixsaN5YCoUHYyn7DmtxUsVp2dHkxhNWMgXqz/Uon6t3Ea8JjSVt4ei5pI6cBFR+0rPzkKSJITD7DJHRK2B72iI2oSqyBgoN2V4VZP9RUdgImU3XI43ErfwwNEJDMezdYEpYqrlrnj1C9dWtsOmtnovjojWjNnZWYTDYcgyG8MQUWtgICIiAOVW4BED/REDaHJZv+MITKZz5UvxSkHpVDk0/eLYFIbj2brFKIO66rYTH4zUtBYvh6aoqS25ExARrT2pVArRaNTrYRARuRiIiOiMyLKE3pCO3pCOSzZFG/YLITCVzlWDUk1wOvBiKTBl89XA5Pcp5Tbi1erSQDk0DUVNdPgZmIjWo9nZWWzctMnrYRARuRiIiGhZSJKE7qCO7qCOi4YaJ0sLITCTyTeswzQ8k8VjJ2Zw95ERZHJF93hTk925S5XQNFQzj6k76GNgIlpjhBBIpVLo6OjweihERK4VC0S33XYb/uu//gtHjhyBz+dDPB5vOObkyZO47rrr8MADDyAYDGLPnj3Yt28fVJU5jWi9kSQJnQEfOgM+XDjYPDAlsoWGDnmlLnkJ/NdTMaSsgnu8rsqloFQzj6m2+UNPUIcsMzARtRLLspDP59HR0en1UIiIXCuWPHK5HN7xjndg165d+Kd/+qeG/cViEW95y1vQ39+Phx9+GKOjo3j3u98NTdPw2c9+dqWGRUQtSpIkRP0aon4N2zc07z6VzOYxkrBwaiZbDk6l0PT0cAL//cwY4tnqatmaIpWaPsxtLV7e7g3pUBiYiFZVKpkEAFaIiKilrFgg+tSnPgUAuOuuu5ru/8lPfoJnn30WP/3pT9HX14eLL74Yn/nMZ3DTTTfhr/7qr+DzLd8icES0PoRNDWFTw/n9oab7Z+0CRuL11aXhuIXnYynsf34cM5lqYFJlCRsiRt0cJrdbXoeJvpAOVWEXLKLllGQgIqIW5Nm1aQcOHMCFF16Ivr4+97arrroK1113HZ555hlccsklXg2NiNaooK7i5X0hvLyveWDK5Ao1LcWroemFiTR+fmwSk7M591hFltAf1usbPtQEp/6IAY2BieisJJJJhMNhaBrb8hNR6/AsEMVisbowBMD9OhaLzXs/27Zh27b7deWvTUREC/H7VLysN4iX9Qab7rfyRYwkrLoOecPxLH4zlcZDL0xhIlX92SNLQF/YwECk3Fo8Wr8e04aIAZ/KwERUK5lIoKu72+thEBHVOatA9PGPfxy33377aY957rnncP755y9pUKezb98+93I8IqLlZGgKzu0O4NzuQNP9dr6I0aSF4RkLI4ksTlWaP8xkcfD4NMZSNkR57VpJAnqCesOleJU5TQMRA4amrOKrI/JeIpFc0fcIRESLcVaB6CMf+Qj+9E//9LTHnHvuuWf0WP39/Th48GDdbWNjY+6++dx8883Yu3ev+3UymcTGjRvP6DmJiJZC1xSc0xXAOV3NA1Ou4GAsWVmw1irPZyo1f3jipThGExYcUT2+J+ir64xXWbR2qMPEQMSE6WNgovWjWCxidjaFru4ur4dCRFTnrAJRT08Penp6luWJd+3ahdtuuw3j4+Po7e0FANx3330Ih8PYvn37vPfTdR26ri/LGIiIlpNPlbGx04+Nnf6m+/NFB2NJu3w5Xv08pieHE4glLBRqElNnQKubwzRUsz0QNRHUuUQBrR3JZBJCCHTzkjkiajEr9tv05MmTmJ6exsmTJ1EsFnHkyBEAwMte9jIEg0G86U1vwvbt23HttdfijjvuQCwWwyc+8Qlcf/31DDxEtC5pioyhDhNDHWbT/UVHYDxlY3gmi1O1i9fGs3h2NInRhIV8sRqYon4NgzVzmOY2fggZnLhOrSM+MwMA6O5enj+sEhEtlxULRJ/85CfxjW98w/260jXugQcewOtf/3ooioJ77rkH1113HXbt2oVAIIA9e/bg05/+9EoNiYiopSnlVuAbIgYuQ2NbYscRmJi1y+sw1Td/eODoBEYSFnIFxz0+bKju3KWBSCksDXVUg1PYUCFJXIuJVkc8HkcwGIRhGF4PhYiozooForvuumveNYgqNm/ejB/96EcrNQQionVFliX0hQ30hQ1c2mS/4whMpXN1azCV5jBl8f9emMJwPAsrXw1MQV11w9FQbYWp3PwhamoMTLRs4vE4esqXyBMRtRJegE5EtE7IsoSekI6ekI6LN0Yb9gshMJ3ONSxcOxLP4pHjUxiOW8jkiu7xfp/idsgbmNNWfDBqoDPgY2CiMzYzM4MLLrjQ62EQETVgICIiahOSJKErqKMrqOOioUjDfiEE4tm820q80vjhVDyLQydm8P1fZpG2q4HJ0GQMRMpBqclaTN1BBiYqsW0b6XTabaJERNRKGIiIiAhAKTB1+H3o8PtwwUC4Yb8QAkmrUA1KlTlMCQu/fCmBHz0VQ9IquMf7VBkDEaPhcryBSKm1eE9QhywzMLWDmXJDhd45C7ITEbUCBiIiIjojkiQhYmqImBpesaExMAFAysq7l+TVdsl7ZjSFnzw3jngm7x6rKaUmErWtxSvzmQY7TPSGdCgMTOvCzPQ0FEVBVxfXICKi1sNAREREyyZkaDi/X8P5/aGm+9N2oa5DXiU8HR1L4f6j45hOVwOTKkvojxgNc5cqi9n2h3WoirxaL42WYHp6Gj09vZBlni8iaj0MREREtGoCuoptvUFs6w023Z/JFTCasOa0Frfw4kQavzg2iYnZnHusIkvoC+nlS/GM0nymcoe8oaiJvrABn8o34K1genoa52zZ4vUwiIiaYiAiIqKW4fep2NoTxNae5oHJzhdLQammQ96peBYnprJ4+IVpTMzaEOW1ayUJbmCqdMkb6jDLwakUoBiYVl4+n0cikUB/f7/XQyEiaoqBiIiI1gxdU7ClO4At3YGm+3MFB6MJq2EtppF4Fo+fmMFY0oIjqsf3hHQMRRvnMVWaPxiaskqvbP2amZmBEAL9/Ru8HgoRUVMMREREtG74VBmbu/zY3OVvuj9fdBBLWhieqQalSmh64qU4YkkbxZrE1B301QelmtbiA1EDfh9/jS5kanISiqKgp6fH66EQETXFn+RERNQ2NEXGxg4/NnY0D0yFooOxlF1ahylR01o8buHp4SRGExYKNYGpM6DVzV2a2/whqPPX7OTkJPr6+qEorLYRUWviT2oiIqIyVZHdUNNM0RGYSNk4Fc/WtxafyWL/aAojiSzyxWpgippaXWe8ytyloXKAChnaar00z0xNTeG88873ehhERPNiICIiIjpDSrkVeH/EwGWbOxr2O47AxKztNnyoncf04LEJjMQt2AXHPT5sqG7Dh7pL8srNHyKmCklau2sxWZaFVCqFgYEBr4dCRDQvBiIiIqJlIssS+sIG+sIGXrkp2rBfCIGpdM5tJ17b/OHhF6YwHM8im68GpoCulEJSxChfllc/n6nDr7V0YJqYmAAADAwOejwSIqL5MRARERGtEkmS0B3U0R3UsXNj434hBGYyeZyqNHyoaf5w8PgMTsVHkMkV3eNNrXqJX6Uz3mCHiaFycOoK+DwNTJMTE/D7/YhEIp6NgYhoIQxERERELUKSJHQGfOgM+HDRYGOIEEIgkS24laWRuOXOZzp8Mo4fzGQxaxfc4w1NxkCkWlUamtNavDvggyyvXGCamJjA0NBQS1exiIgYiIiIiNYISZIQ9WuI+jXsGAg3PSaZzdd1yKsEp6eGE/jxMzEkstXA5FPlUlWppjNeZVvL5+rWbDpbxWIRk5OT2L7jgsU/CBHRKmAgIiIiWkfCpoawqeH8/lDT/bNWAcOJrNtafLgcnJ4dTeG+58Yxk8m7x0p4Jb5/zwQ6/XFEDAVRU637HDFURE3FvS2sK1DKFafp6WkUi0VsHBpalddNRLRYDERERERtJGioOM8I4by+5oEpbRcwmrDw9PERfO8nD0LWI7CLGjJJGcdnZGSK1Q8Hcy+FEzBlAb/iwJDykIvb8OL9MUjS2Mq/MCLyhOQUcI7Xg1giBiIiIiJyBXQVL+sNYkvXVgRTJ5BIJAAUyx9VQgCWA6QLMtIFyf3IFCSkixLSBQVCiUJTuSAr0XomSbLXQ1gyBiIiIiJqoCgKdu/+ba+HQUQtTgiBfP4iAICmrc3FphmIiIiIiIhoUSRJgs/n83oYS7L2a1xERERERESLxEBERERERERti4GIiIiIiIja1pqfQyREadW4VCrl8UiIiIiIiMhLlUxQyQhnYs0HosqL3nruFo9HQkRERERErSCVSiESiZzRsZI4m/jUghzHwcjICEKhECRp7gJxdLaSySQ2btyIl156CeFw2OvhtCWeA+/xHHiP56A18Dx4j+fAezwH3jubcyCEQCqVwsDAAGT5zGYHrfkKkSzLGBoa8noY6044HOZ/eo/xHHiP58B7PAetgefBezwH3uM58N6ZnoMzrQxVsKkCERERERG1LQYiIiIiIiJqWwxEVEfXddx6663Qdd3robQtngPv8Rx4j+egNfA8eI/nwHs8B95b6XOw5psqEBERERERLRYrRERERERE1LYYiIiIiIiIqG0xEBERERERUdtiICIiIiIiorbFQESur3zlKzjnnHNgGAauuOIKHDx40OshrVv79u3Dq171KoRCIfT29uLtb387jh49WneMZVm4/vrr0dXVhWAwiD/8wz/E2NiYRyNe/z73uc9BkiTceOON7m08B6tjeHgYf/Inf4Kuri6YpokLL7wQjz/+uLtfCIFPfvKT2LBhA0zTxO7du3Hs2DEPR7y+FItF3HLLLdiyZQtM08TWrVvxmc98BrU9l3gOltfPf/5z/P7v/z4GBgYgSRK+973v1e0/k+/39PQ0rrnmGoTDYUSjUbzvfe/D7OzsKr6Kte105yCfz+Omm27ChRdeiEAggIGBAbz73e/GyMhI3WPwHCzNQv8Pan3gAx+AJEn44he/WHf7cp0DBiICAHz729/G3r17ceutt+Lw4cPYuXMnrrrqKoyPj3s9tHXpwQcfxPXXX49HHnkE9913H/L5PN70pjchnU67x3z4wx/GD3/4Q3znO9/Bgw8+iJGREVx99dUejnr9euyxx/AP//APuOiii+pu5zlYeTMzM7jyyiuhaRruvfdePPvss/ibv/kbdHR0uMfccccd+NKXvoSvfvWrePTRRxEIBHDVVVfBsiwPR75+3H777bjzzjvx93//93juuedw++2344477sCXv/xl9xieg+WVTqexc+dOfOUrX2m6/0y+39dccw2eeeYZ3Hfffbjnnnvw85//HO9///tX6yWseac7B5lMBocPH8Ytt9yCw4cP47vf/S6OHj2Kt771rXXH8RwszUL/DyruvvtuPPLIIxgYGGjYt2znQBAJIS6//HJx/fXXu18Xi0UxMDAg9u3b5+Go2sf4+LgAIB588EEhhBDxeFxomia+853vuMc899xzAoA4cOCAV8Ncl1KplNi2bZu47777xOte9zrxoQ99SAjBc7BabrrpJvHa17523v2O44j+/n7x+c9/3r0tHo8LXdfFv//7v6/GENe9t7zlLeK9731v3W1XX321uOaaa4QQPAcrDYC4++673a/P5Pv97LPPCgDisccec4+59957hSRJYnh4eNXGvl7MPQfNHDx4UAAQJ06cEELwHCy3+c7BqVOnxODgoHj66afF5s2bxd/+7d+6+5bzHLBCRMjlcjh06BB2797t3ibLMnbv3o0DBw54OLL2kUgkAACdnZ0AgEOHDiGfz9edk/PPPx+bNm3iOVlm119/Pd7ylrfUfa8BnoPV8oMf/ACXXXYZ3vGOd6C3txeXXHIJvva1r7n7jx8/jlgsVnceIpEIrrjiCp6HZfKa17wG+/fvx69+9SsAwC9/+Us89NBDePOb3wyA52C1ncn3+8CBA4hGo7jsssvcY3bv3g1ZlvHoo4+u+pjbQSKRgCRJiEajAHgOVoPjOLj22mvx0Y9+FDt27GjYv5znQF3yaGnNm5ycRLFYRF9fX93tfX19eP755z0aVftwHAc33ngjrrzySlxwwQUAgFgsBp/P5/7grejr60MsFvNglOvTt771LRw+fBiPPfZYwz6eg9Xx4osv4s4778TevXvxF3/xF3jsscfwwQ9+ED6fD3v27HG/181+PvE8LI+Pf/zjSCaTOP/886EoCorFIm677TZcc801AMBzsMrO5Psdi8XQ29tbt19VVXR2dvKcrADLsnDTTTfhXe96F8LhMACeg9Vw++23Q1VVfPCDH2y6fznPAQMRkceuv/56PP3003jooYe8Hkpbeemll/ChD30I9913HwzD8Ho4bctxHFx22WX47Gc/CwC45JJL8PTTT+OrX/0q9uzZ4/Ho2sN//Md/4Jvf/Cb+7d/+DTt27MCRI0dw4403YmBggOeA2l4+n8cf/dEfQQiBO++80+vhtI1Dhw7h7/7u73D48GFIkrTiz8dL5gjd3d1QFKWhe9bY2Bj6+/s9GlV7uOGGG3DPPffggQcewNDQkHt7f38/crkc4vF43fE8J8vn0KFDGB8fxytf+UqoqgpVVfHggw/iS1/6ElRVRV9fH8/BKtiwYQO2b99ed9srXvEKnDx5EgDc7zV/Pq2cj370o/j4xz+Od77znbjwwgtx7bXX4sMf/jD27dsHgOdgtZ3J97u/v7+h6VGhUMD09DTPyTKqhKETJ07gvvvuc6tDAM/BSvvFL36B8fFxbNq0yf0dfeLECXzkIx/BOeecA2B5zwEDEcHn8+HSSy/F/v373dscx8H+/fuxa9cuD0e2fgkhcMMNN+Duu+/G/fffjy1bttTtv/TSS6FpWt05OXr0KE6ePMlzskze+MY34qmnnsKRI0fcj8suuwzXXHONu81zsPKuvPLKhpbzv/rVr7B582YAwJYtW9Df3193HpLJJB599FGeh2WSyWQgy/VvBxRFgeM4AHgOVtuZfL937dqFeDyOQ4cOucfcf//9cBwHV1xxxaqPeT2qhKFjx47hpz/9Kbq6uur28xysrGuvvRZPPvlk3e/ogYEBfPSjH8V///d/A1jmc7C4XhC03nzrW98Suq6Lu+66Szz77LPi/e9/v4hGoyIWi3k9tHXpuuuuE5FIRPzsZz8To6Oj7kcmk3GP+cAHPiA2bdok7r//fvH444+LXbt2iV27dnk46vWvtsucEDwHq+HgwYNCVVVx2223iWPHjolvfvObwu/3i3/91391j/nc5z4notGo+P73vy+efPJJ8ba3vU1s2bJFZLNZD0e+fuzZs0cMDg6Ke+65Rxw/flx897vfFd3d3eJjH/uYewzPwfJKpVLiiSeeEE888YQAIL7whS+IJ554wu1gdibf79/5nd8Rl1xyiXj00UfFQw89JLZt2ybe9a53efWS1pzTnYNcLife+ta3iqGhIXHkyJG639O2bbuPwXOwNAv9P5hrbpc5IZbvHDAQkevLX/6y2LRpk/D5fOLyyy8XjzzyiNdDWrcANP3453/+Z/eYbDYr/vzP/1x0dHQIv98v/uAP/kCMjo56N+g2MDcQ8Rysjh/+8IfiggsuELqui/PPP1/84z/+Y91+x3HELbfcIvr6+oSu6+KNb3yjOHr0qEejXX+SyaT40Ic+JDZt2iQMwxDnnnuu+Mu//Mu6N348B8vrgQceaPo7YM+ePUKIM/t+T01NiXe9610iGAyKcDgs3vOe94hUKuXBq1mbTncOjh8/Pu/v6QceeMB9DJ6DpVno/8FczQLRcp0DSYiapaiJiIiIiIjaCOcQERERERFR22IgIiIiIiKitsVAREREREREbYuBiIiIiIiI2hYDERERERERtS0GIiIiIiIialsMRERERERE1LYYiIiIiIiIqG0xEBERERERUdtiICIiIiIiorbFQERERERERG2LgYiIiIiIiNrW/w83VTWq1A303wAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 1000x400 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "lens.draw(num_rays=5)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
